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Phase variation, or stochastic switching between alternative states of gene expression, is common among mi- 
crobes, and may be important in coping with changing environments. We use a theoretical model to assess 
whether such switching is a good strategy for growth in environments with occasional catastrophic events. We 
find that switching can be advantageous, but only when the environment is responsive to the microbial popula- 
tion. In our model, microbes switch randomly between two phenotypic states, with different growth rates. The 
environment undergoes sudden "catastrophes", the probability of which depends on the composition of the pop- 
ulation. We derive a simple analytical result for the population growth rate. For a responsive environment, two 
alternative strategies emerge. In the "no switching" strategy, the population maximises its instantaneous growth 
rate, regardless of catastrophes. In the "switching" strategy, the microbial switching rate is tuned to minimise the 
environmental response. Which of these strategies is most favourable depends on the parameters of the model. 
Previous studies have shown that microbial switching can be favourable when the environment changes in an 
unresponsive fashion between several states. Here, we demonstrate an alternative role for phase variation in 
allowing microbes to maximise their growth in catastrophic responsive environments. 
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I. INTRODUCTION 

Microbial cells often exhibit reversible stochastic switching 
between alternative phenotypic states, resulting in a heteroge- 
neous population. This is known as phase variation |1-3|. 
A variety of molecular mechanisms can lead to phase varia- 
tion, including DNA inversion, DNA methylation and slipped 
strand mispairing [T, T|. These are generally two-state sys- 
tems without any underlying multistability |4, 5|; however, 
bistable genetic regulatory networks can also lead to stochas- 
tic phenotypic switching |6-9|. The biological function of 
phase variation remains unclear, but it has been suggested that 
it can allow microbes to evade host immune responses, or ac- 
cess a wider range of host cell receptors |3, 10|. Theoretical 
work has focused on phase variation as a mechanism for cop- 
ing with environmental changes. According to this hypothe- 
sis, a fraction of the population is maintained in a state which 
is currently less favourable, but which acts as an "insurance 
policy" against future environmental changes [1 Ij. 

In this paper, we present a theoretical model for switching 
cells growing in an environment which occasionally makes 
sudden attacks on the microbial population. Viewing the sit- 
uation from the perspective of the microbes, we term these 
catastroplies. These catastrophes affect only one phenotypic 
state. Importantly, the environment is responsive: the catas- 
trophe rate depends on the microbial population. By solving 
the model analytically, we find that there are two favoured 
tactics for microbial populations in environments with a given 
"feedback function": keep all the population in the fast grow- 
ing state, regardless of the environmental response, or alter- 
natively, use switching to maintain a population balance that 



reduces the likelihood of an environmental response. Which 
of these strategies is optimal depends on the parameters of 
the model. In the absence of any feedback between the 
population and environment, phase variation is always un- 
favourable. However, as the environment becomes more re- 
sponsive switching can be advantageous. 

Previous theoretical studies have considered models in 
which the environment flips randomly or periodically between 
several different states, each favouring a particular phenotype. 
The case of two environmental states and two phenotypes has 
been well-studied lfT2l - [r7l . This work has shown that the to- 
tal growth rate of the population can be enhanced by phe- 
notypic switching (compared to no switching), for some pa- 
rameter regimes, and that the optimum switching rate is tuned 
to the environmental flipping rate. Several studies have also 
compared random switching to a strategy where cells detect 
and respond to environmental changes. Wolf et al. {\J \ used 
simulations to show that in this case the advantage of random 
switching depends on the accuracy of environmental sensing, 
while in a theoretical study Kussell and Leibler [18 | showed 
that the advantages of random switching depend on the cost of 
environmental sensing for a model with n phenotypic states 
and n different environments. The predictions of the "two en- 
vironment, two phenotypic state" model have recently been 
verified experimentally with a tunable genetic switch in the 
yeast Saccharomyces cerevisiae |19|. 

Here, we consider a different scenario to the above- 
mentioned body of work. Rather than considering multiple 
environmental states, our model has a single environment, 
which undergoes occasional, sudden and instantaneous catas- 
trophes. We assume that the more slowly-growing micro- 
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bial phenotypic state is resistant to these catastrophes. Catas- 
trophic events are Hkely to be a common feature of microbial 
population dynamics in nature. For example, microbes infect- 
ing an animal host may be subject to sudden "flushing" due to 
diarrhoea or urination, to which they may be resistant if they 
are able to attach to the wall of the host's intestinal or urinary 
tract. Another example of a catastrophe might be sudden ex- 
posure of a population to antibiotics: here, cells that are in the 
non-growing "persister" state survive, while others are killed 
II20I |2i I . We do not, however, aim to model a specific bio- 
logical case, but rather to construct a generic model leading to 
general conclusions. 

Importantly, and in contrast to previous models, we include 
in our model feedback between the microbial population and 
the environment: the probability of a catastrophe depends on 
the state of the population. Although our model is very gen- 
eral, many examples exist in nature in which environmental 
responses are triggered by characteristics of a growing micro- 
bial population, the most obvious perhaps being a host im- 
mune response |22|. Our work leads us to propose an alter- 
native possible role for phase variation, to our knowledge not 
considered in previous theoretical work: we find that in re- 
sponsive catastrophic environments, switching can allow the 
population to maximise its growth rate while minimising the 
environmental response. 

The paper is organised as follows. In sectionjll] we present 
our model. We derive an analytical result for the steady-state 



statistics of the model in Section III and we use this to predict 
the optimal strategies for microbial growth, as a function of 



the model parameters, in Section IV Finally, we present our 
conclusions in Section IVl 



II. MODEL 

We consider two microbial sub-populations A and B, rep- 
resenting two different phenotypic states. Between catastro- 
phes, microbes in these sub-populations grow exponentially at 
rates ja and jb, and switch between states with rates kA and 
fcs (A to B and _B to ^ respectively). However, this growing 
regime can be ended suddenly by a catastrophe, which con- 
sists of a sharp decrease in the size of the A sub-population. 
After the catastrophe, the population dynamics restarts. 

Between catastrophes the dynamics of the numbers of mi- 
crobes UA and in the two sub-populations are defined by 
the following system of differential equations 



duA 



dns 
~dt 



"iAnA + ksriB - kAUA 



iBnB + kAflA - kBTlB 



(la) 



(lb) 



This description assumes that the population sizes ua and ub 
are large enough to be considered as continuous variables. We 
assume that 7^ > 75, which means that the A sub-population 
proliferates faster than the B sub-population. 

Whenever a catastrophe takes place, the population size ua 
drops instantaneously to some new value < ua, with a 



probability i/(n^|n^). The rate at which catastrophes happen 
depends on the population size through an environmental re- 
sponse function f3{nA,nB)- This function characterises the 
rate at which the environment responds to the growing popu- 
lation. The two functions /3 and i' are discussed in detail in 
Section [ILB] A typical trajectory for the sizes of the A and B 
sub-populations, for a particular choice of /3 and i^, is shown 
in the top panel of Fig{T] 



A. Fitness 

As shown in ||T4l, the two-variable system defined by Eqs. 
([T]l can be replaced by a non-linear dynamical equation for a 
single variable. This variable, /, is the fraction of the total 
population in the A state: 



fit) 



nA 

nA + nB 



(2) 



If we consider the dynamics of the total population n{t) 
UAit) + nB{t), then, from ([TJ it follows that Ifl4l : 



dn{t) 
dt 



lAnA + IsnB = {iB + Mf)n{t) 



(3) 



where A7 — ^a — 1b > 0. The above equation shows that 
/ is linearly related to the instantaneous growth rate of the 
population (which is given by 7^ + A7/). For this reason, 
and following [ 14,1 , we refer to / as the population fitness. 

The dynamical equation for the population fitness can be 
determined from Eqs.([T|i and corresponds to: 



dt 



(4) 



where we define v{f) as the time evolution function for the 
fitness, and f± are the two roots of the quadratic equation: 
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A7 
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A7 







(5) 



One can check that the smaller root takes values /_ < 0, 
while the larger root takes values < /+ < 1. Hence, the 
population fitness increases towards a plateau value /+, until 
a catastrophe happens, upon which it is reset to a lower value. 
A typical time trajectory for the population fitness is plotted 
in the bottom panel of Fig[T] The time evolution of / is de- 
terministic except at some specific time points (catastrophes) 
where it undergoes random jumps. This model can therefore 
be considered to be a Piecewise Deterministic Markov Pro- 
cess 123] El. 



B. Catastrophes 

The catastrophes in our model have two characteristics: the 
rate at which they happen and their strength {i.e. how many 
microbes are killed). The rate at which catastrophes occur, 
or their probability per unit time, is defined by a feedback 
function /?(/), which we take to depend only on the fitness 
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FIG. 1 Typical trajectory of the system. Top: time evolution 
of the populations ua and ub (semi-log scale). Bottom: Corre- 
sponding time evolution of the fitness f(t) (the fraction of cells 
in the A state). Parameters are: (7^, 7s, ^a, fcs, /*, /^o, ct) = 
(6,1,1,1,0.25,1,-0.75). 



of the population and not on the absolute population size (we 
shall return to this assumption in the Discussion). The func- 
tion /3(/) characterises the response of the environment to the 
growth of the population. When (5 — 0, there are no catastro- 
phes and the fitness will reach the plateau value /+ and stay 
there for ever Nonzero constant values of /? correspond to a 
non-responsive environment in which the catastrophes follow 
Poisson statistics. We shall consider the case of a responsive 
environment characterised by a response function /?(/) which 
depends on the population fitness. In particular, we consider 
a nonlinear response function that has a sigmoid shape. Thus, 
the probability per unit time of a catastrophe is very low when 
the population fitness is low, but increases significantly if the 
fitness exceeds some threshold value. This scenario might cor- 
respond to a detection threshold in the environment's sensitiv- 
ity to population growth. 

The precise environmental response function that we con- 
sider is the following: 



/3a(/) 




i-r 



+ {f-.f*y 



(6) 



Although this function is defined for all —00 < / < 00, the 
relevant interval for the fitness is < / < 1. Typical shapes 
for this function are shown in Fig|2] The parameter ^ is the 
asymptotic value of /?> when / is large, and we refer to ^ as 
the saturated catastrophe rate. As the population fitness / in- 
creases, /3\ increases from to ^ around the threshold value 
/* at which = ^/2. Finally, the parameter A determines 
the sharpness of the threshold. For small values of A the func- 
tion /3x (/) approaches a step function 

M) = mf-r)- (7) 

As A increases, the function broadens and becomes linear over 




FIG. 2 The response function Px (/), plotted for different values of 
A, with 5 = land/* = 1/2. 



a range of / near /* : 



if-n 

A 



O 



(8) 



while when the parameter A becomes very large (A — 00), 
Pxif) becomes constant (independent of /) so that the catas- 
trophes become a standard Poisson process with parameter 
e/2: 



/3oo(/) = e/2. 



(9) 



We emphasise that we have chosen this particular sigmoid 
function (j6]) as the A parameter allows a convenient tuning of 
its shape and thus the degree of environmental responsiveness. 
However, our conclusions are not affected by the particular 
choice of sigmoid function. 

We now turn to the function describing the catastrophe 
strength, i/(n^|n^). This is the probability that, given that 
UA cells of type A are present before the catastrophe, will 
remain after the catastrophe. In order to retain our descrip- 
tion of the model in terms of the population fitness, we shall 
consider that v only depends on only through the ratio 
n'j^/nA- Then the normalisation of v implies that 



v{nA\nA) — TijiAluA) 

UA 



(10) 



where dxT[x) — 1. When a catastrophe happens, the pop- 
ulation size is reduced by a random factor sampled from the 
distribution [i.e. the new size = x u where ua 
is the size before the catastrophe, and u is a random number 
(0 < M < 1) sampled from the distribution J^]. This allows us 
to associate to each jump — > a fitness jump / — > /', 
where /' = n^/(n^ + ub)- The size of these jumps will be 
distributed according to: 



1-/ 



/'(!-/) 

f{i-nj [i-f'ff 



(11) 



Eq.dli can be obtained by rewriting expression ( [T0| for 
v{n'j^\nA) as a function of / and /' and including the Jacobian 
of the transformation. 

In this paper, we shall consider the simple case where 
J-{x) = (a + l)x°', with a > — 1. The explicit expression for 
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FIG. 3 The jump distribution v(rij^\nA), plotted for various values 
of Q. 



v(r)!\n) thus reads: 

, (a + 1) /^n'^" 
v(n \n) — ' ' 



a > -1 



(12) 



This choice is made primarily in order to allow us to solve 
the model analytically: the choice implies that fJ.{f'\f) fac- 
torises (see ( 13 1) which then allows the integral equation for 
the probability flux balance ([TSj to be solved. Moreover, the 
choice of a power law distribution for (/(n^jn^) is general in 
that allows for increasing, decreasing or flat functional forms. 
The function is plotted in Figj3]for various values 

of a. For negative a values, the distribution is biased towards 
far-reaching catastrophes that reduce fitness significantly. The 
case a = corresponds to jumps sampled from a uniform 
distribution, whereas positive values give a distribution biased 
towards weaker catastrophes. The parameter a can therefore 
be used to tune the strength of the catastrophes (although in 
this work we shall always consider negative a values, corre- 
sponding to strong catastrophes). We note here that with our 
choice for J-{x) the jump distribution can be expressed as 



df m{f) 

where m(0) = 0, / df^i{f\f) ^ 1, and with, 

/ 



(13) 



m{f) 



1-/ 



(14) 



III. STEADY-STATE STATISTICS 

We now derive the steady-state probability distribution for 
the population fitness, p{f). The distribution p{f) must sat- 
isfy a condition of balance for the probability flux. This con- 
dition reads: 



vimj) 



f 



f 



df dr(3{f)p{funf) . (15) 



The left-hand side of the above equation corresponds to the 
deterministic probability flux due to population growth as de- 
fined in Eq.(|4|) [f{t) increases in time as the population grows. 
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FIG. 4 Illustration of the flux balance condition. The positive proba- 
bility flux due to population growth must be balanced by the negative 
flux due to catastrophes. 



as shown in Figure[T]. The right-hand side describes the prob- 
ability flux arising from catastrophes. In this model, catas- 
trophes always reduce the population fitness. The probabil- 
ity flux due to catastrophes therefore contains contributions 
from all possible jumps that start at some /' > / and end 
at some /" < /. These contributions must be weighted by 
/3(/')p(/'): the probability of having fitness /' and under- 
going a catastrophe. This balance between the fluxes due to 
growth and catastrophes is illustrated schematically in Figure 

El 

Inserting the form ( 13 1 of fi{f'\f) in Eq.( 15 i, the zero flux 
condition becomes: 



(16) 



We now divide the above equation by r7i(/) and take the first 
derivative with respect to /. This yields, in terms of the func- 
tion G = vp/m: 



dG 



f3G 

V 



(17) 



The above differential equation is then easily solved for G. 
The result forp(/), using the form Eq.( 14 1 of m( f) , is finally 



Pif) 



C 

W) 



f 

1-/ 



exp 



(18) 



where C is a normalisation constant. The integral in Eq.( 18 1 
can be performed analytically for the model defined in the 
previous section. The result, which is rather cumbersome, is 
given in the Appendix. 

We present in Fig|5](top panel) some resulting shapes for 
the probability distribution p{f) in the case A = 0, corre- 
sponding to a step function for the environmental response. 
We consider two different values of A7, in each case for fc^ = 
(no switching) and a non-zero switching rate kj^ — kj^* de- 
fined such that /+ = /* (see section IV 1. In these plots we 
see that singularities in p{f) can arise at / = 0, /* or /+ in 
different cases. In particular, the kA — case produces a cusp 
at / = /* (due to the singular nature of the step function /3(/) 
at /*) for small A7, and/or a divergence at / = /+ for large 
A7. On the other hand fc^ = fc^* produces a divergence at 
/+ = /* in both cases. 

Fig|5] (bottom panel) plots trajectories of the fitness corre- 
sponding to the parameter values of Fig|5] These trajectories 
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FIG. 6 Average population fitness (/) as a function of the A ^ B 
switching rate kA, for several values of parameter A. For large values 
of A the average fitness peaks for zero switching, fc^ — 0. However, 
as A decreases, a second peak arise at non zero kA, indicating that 
random switching can be a favourable strategy when the environment 
is responsive. Parameters values are: fcs = 0.1, A7 = 1, a = 
-0.99, 5 = 1 and /* = 0.5. 



FIG. 5 Top panel: Examples of steady-state fitness distribution p{f) 
for parameter values A7= 0.1 (left) and 100 (right); other parame- 
ters are A = 0, fcs = 0.5, = 1, a = -0.99 and /* = 0.75. In 
each plot, the solid line corresponds to a value of fc^ = (no switch- 
ing), while the dashed plot is for kA = fc,4 (switching rate given by 
Equation 1 19[l. Bottom panel: Examples of fitness trajectories corre- 
sponding to the parameter values of the top panel. 



reveal the interplay between two timescales: the time to re- 
lax to the plateau value /* in the absence of catastophes and 
the typical time between catastrophes. The former decreases 
with A7 and the latter is given by 1/^ where ^ is the plateau 
value of the response function f3. A divergence of p{f) at 
f = f+ arises when the plateau value is typically reached 
before a catastrophe occurs. 



IV. OPTIMAL STRATEGIES: TO SWITCH OR NOT TO 
SWITCH? 

The key question to be addressed in this work is whether 
random switching is advantageous to the microbial popula- 
tion in our model. To answer this question, we take advantage 
of the analytical solution 



18 1 to investigate how the time- 



averaged population fitness depends on the rate kA of switch- 
ing from the fast-growing state A to the slow-growing state B. 
We are particularly interested in the effect of the parameter A, 
which controls the sharpness of the environment's response to 
the population. 

In Fig|6]we plot the average population fitness against fc^ 
for several values of A. For a non responsive environment [i.e. 
in the limit of large A, where the catastrophe rate takes the con- 
stant non-zero value (j9])], the population fitness has only 
one (boundary) maximum for switching rate 0. This 

means that the optimal rate of population growth is achieved 
when the bacteria do not switch away from the fittest state 
A. It should be noted that we consider the limit kA — >■ 0, so 
that the population always contains some small residual frac- 
tion in the unfit B state, which becomes a finite fraction of 



the population after a catastrophe. Subsequently, in between 
the catastrophes, the A subpopulation grows quickly to dom- 
inate the population and the fitness evolves towards the value 
/+ = 1, which follows from Eq.(|5]l when kA — 0. 

In contrast, as the environment is made responsive by de- 
creasing the parameter A, a local maximum appears in the 
population fitness, for nonzero switching rate fc^. This im- 
plies that for responsive environments, switching into the 
slow-growing state represents an optimal strategy for the mi- 
crobes. The height of the peak at fc^ 7^ can surpass that of 
the peak at fc^ = 0, showing that random switching can be ad- 
vantageous compared to keeping the whole population in the 
fast-growing state, if the environment is responsive. Thus the 
two maxima coiTespond to two alternative strategies which we 
term switching for the peak at kA = kA and non-switching for 
the peak at kA = 0. 

To gain further insight into the meaning of these two strate- 
gies, and to determine which circumstances favour one strat- 
egy over the other, we focus on the limiting case A = 0, where 
the response function is a step function with its threshold at 
f — f*. We assume that the environmental threshold /* 
is less than the maximum population fitness /+. (If this is 
not the case, the unrealistic situation arises where the popula- 
tion never has a high enough fraction of A cells to trigger any 
catastrophes.) Since depends on the switching rate kA via 
Eq.dSll, this condition /* < implies a maximum value fc^ 
for ka- 



h* — 

Ka — 



{1 - n{Aj f* + kB) 



(19) 



FigjTjshows two examples of how the average fitness (/) de- 
pends on kA in the range Q io k\. One can see that there 
are always two boundary maxima located at A:^ = and at 
kA = k*A\ these correspond to the non-switching and switch- 
ing strategies. 

We plot in Fig[8]typical trajectories of the population fitness 
for the two cases corresponding to the black circles in Fig|7] 
(left panel). These trajectories have the same time-averaged 
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FIG. 7 Average population fitness {/) as a function of kA for a 
threshold response function. The left plot is an instance where the 
maximal fitness occurs at non-zero kA, while the right plot shows an 
example where the maximum fitness occurs at kA ~ 0. In both plots 
parameters are a = —0.99, /* = 0.75, Po — 1 and fcs — 0.8; 
A7 = 1 (left) and A7 = 4 (right). 



FIG. 8 Two typical trajectories illustrating the two possible strate- 
gies: in the left panel, the switching rate kA = 0, the fitness evolves 
towards its maximal value / = 1 but many catastrophes are in- 
duced; in the right panel fc^ is set to a value smaller but very close 
to k\ {kA ~ 0.6155 — 0.997fcJ), the fitness evolves to just above 
the threshold value with relatively few catastrophes induced. These 
strategies correspond to the two black points in the top left plot of 
Fig|7] The value of /* is 0.75, as indicated by the red dashed line, 
and the average fitness is the same in both plots ((/) = 0.652585). 



population fitness, but they show very different dynamical be- 
haviour The non-switching strategy (fc^ = 0) is characterised 
by a fast evolution of the fitness towards its maximum = 1. 
However this triggers frequent catastrophes which cause sud- 
den decreases in fitness. In contrast, for the switching strat- 
egy (kA k\), the fitness has a slower growth towards a 
plateau value at the detection threshold /*. In this way, the 
population reduces the frequency of catastrophes by maintain- 
ing itself in a heterogeneous state with a non-zero fraction of 
slower-growing cells which do not trigger catastrophes. 

We next consider how the parameters of our model affect 
the balance between the switching and non-switching strate- 
gies. To this end, we plot "phase diagrams" showing which of 
these two strategies achieves a higher population growth rate 
for a given set of parameters. Figure |9] considers the param- 
eters describing the microbial population: the difference A7 
in growth rate between the A and B states and the switching 
rate fc^ from the slow-growing B state to the fast-growing A 
state. This diagram shows that the switching strategy is only 
favourable when the B state does not caiTy too high a cost in 
terms of growth rate ( A7 not too large) and when switching to 
the B state is unlikely to be immediately followed by a reverse 
switch back into the A state (kg not too large). In Fig 10 we 
consider instead the parameters describing the environmental 
response: the detection threshold /* and the saturated catas- 
trophe rate ^. Here, we see that the switching strategy (i.e. 
attempting to avoid catastrophes) is favoured when the satu- 
rated catastrophe rate ^ is high or when the threshold value 
/* is high (since for high thresholds the population does not 
have to pay a very high price in terms of B cells to avoid 
triggering catastrophes). For very low detection thresholds 
/*, lower than typical values of the fitness, the environmen- 
tal response will almost always "detect" the population, and 
the environmental behaviour will thus be similar to the situa- 
tion of a non-responsive environment, which corresponds to 
the limiting case where /* = 0. In this case, as discussed 



1.5 



cq 1-0 

a; 



0.5 
0.0 









non-switching 




favoured 


switching 




favoured 


1 1 1 


) 1 2 


345 




A7 



earlier, non-switching is the optimal strategy. Figure 10 also 
demonstrates the effect of changing the catastrophe strength 
parameter a (dashed and dotted lines). The switching strat- 
egy is favoured by strong catastrophes (negative a), while the 



FIG. 9 "Phase diagram" showing the optimal strategy as a function 
of the parameters fcs and A7 which describe the microbial popula- 
tion. In the shaded area switching is the optimal strategy, while in 
the unshaded area, the non-switching strategy leads to a larger av- 
erage population fitness. The other parameters are: a — —0.99, 
r =0.75,C = 1. 



non-switching strategy is more likely to be optimal for weak 
catastrophes (i.e. larger positive a). All this points to the 
conclusion that in general, switching tends to be an advan- 
tageous strategy when the characteristics of the catastrophic 
environment are particularly adverse (large ^ and negative a) 
and when the detection threshold is not too low. 



V. DISCUSSION AND FURTHER DIRECTIONS 

In this work we have considered the possible advantages of 
phase variation (random switching between phenotypic states) 
for a microbial population in a catastrophic responsive envi- 
ronment. To this end, we solved analytically for the steady- 
state statistics of a model which includes two microbial sub- 
populations that grow and switch, and a single environment 
which occasionally mounts catastrophic attacks on the micro- 
bial population. Importantly, the model includes feedback be- 
tween the state of the population and the frequency of catas- 
trophic events via an environment response function which 
depends on the population through its fitness i.e. the instanta- 
neous rate of growth. Our results show that, when the environ- 
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FIG. 10 "Phase diagram" showing the optimal strategy as a func- 
tion of the parameters ^ and /* which describe the environmen- 
tal response. In the shaded area switching is the optimal strategy, 
while in the unshaded area, the non-switching strategy leads to a 
larger average population fitness. The other parameters are A7 = 1 
and fes = 1. The black solid line shows the "phase boundary" 
for a = —0.5, while the dotted and dashed line correspond to 
a = -0.99 and q = 0. 



ment is responsive to the population, switching can increase 
the average fitness (i.e. growth rate) of the population. A gen- 
eral picture emerges from our work of two competing strate- 
gies for dealing with a catastrophic responsive environment: 
not switching and thus maximising the instantaneous growth 
rate regardless of catastrophes versus using switching to tune 
the population to reduce the likelihood of catastrophes. 

An important feature of this work is the fact that we are able 
to solve the model analytically, leading to an explicit formula 
for the population fitness as a function of the model parame- 
ters. To achieve this analytical result, we make a number of 
assumptions, the most important being that the environmental 
feedback depends on the instantaneous growth rate rather than 
on the population size. Although this is a rather idealized as- 
sumption, microbe-host interactions are in reality likely to be 
sensitive to microbial growth rate |25|, since several intracel- 
lular small molecules and proteins, including ppGpp, cAMP 
and H-NS, whose concentrations are growth-rate dependent 
|l26 27 1, have been shown to regulate microbial virulence fac- 
tors I28H31. 

The main conclusion of our work is that phase variation 
can provide a mechanism by which a microbial population 
can tune its composition so as to minimize the likely environ- 
mental response, thus increasing its average growth rate (or 
average fitness). The model thus provides an alternative sce- 
nario for the role of phase variation to those proposed in other 
theoretical studies, which we now take the opportunity to re- 
view briefly. 

Various works have considered models in which the envi- 
ronment flips randomly or periodically between several differ- 
ent states, each favouring a particular cell phenotype. These 
models do not include feedback between the population and 
the environmental flipping rate. For the case of two envi- 
ronmental states and two cellular phenotypes, Lachmann and 
Jablonka |12| considered a discrete time model with a pe- 
riodic environment, while Ishii et al ifTSll addressed a sim- 
ilar problem but explicitly looked for the evolutionary sta- 
ble state. Thattai and Van Oudenaarden |14| also considered 



the two-environment, two-phenotype case, using a continu- 
ous time model with Poissonian switching of the environ- 
ment, while a detailed analytical treatment of this case was 
presented by Gander et al. iTTSl and a simulation study was 
carried out by Ribeiro lfT6l with a more detailed model of 
the phenotypic switching mechanism, while Wolf et al. ifTTll 
simulated a model that also included environmental sensing. 
These studies showed that the total growth rate of the popula- 
tion can be enhanced by phenotypic switching (compared to 
no switching), for some parameter regimes, and that the op- 
timum switching rate is tuned to the environmental flipping 
rate. A similar model, but aimed specifically at the case of the 
persister phenotype, in which cells grow very slowly but are 
resistant to antibiotics ll20l . was considered by Kussell et al. 
for a periodic environment f2T|. In this model, the growth rate 
of the non-persister phenotype is negative (signifying popula- 
tion decrease) in the "antibiotic" environment. 

Several other studies have considered random switching 
from a different context: as a means to avoid the need for sens- 
ing and responding to environmental changes, in the case that 
environmental sensing is inaccurate, faulty or expensive. In 
this context, Kussell and Leibler 1 18| considered theoretically 
a model with many environments and many cellular states, 
where a cost is attached to sensing environmental changes, 
while Wolf et al. |T7| simulated a two state-two environment 
model where sensing was subject to a variety of possible de- 
fects. Both these studies concluded that random switching can 
be a good strategy to overcome disadvantages associated with 
environmental sensing. 

In a somewhat different approach. Wolf et al. 131] used 
simulations to study a two state-two environment model in 
which the growth rate of the A and B states is frequency 
dependent — i.e. a given microbial sub-population grows 
faster when its abundance is low. Such "frequency-dependent 
selection" is well known to promote population heterogene- 
ity; however. Wolf et al. did not find any advantages for re- 
versible switching as a means to generate this heterogeneity 
as opposed to terminal cellular differentiation. In a sense, the 
model presented in this paper also incorporates frequency de- 
pendent selection, since catastrophes are less likely when the 
A sub-population is small. However, in contrast to Wolf et al., 
we find that reversible switching does play an important role. 
If switching in our model were not reversible, there would be 
no way for the fast-growing A subpopulation to regenerate 
from the surviving B cells after a catastrophe. 

Although the majority of theoretical work in this area, in- 
cluding that presented in this paper, has focused on the inter- 
play between cellular switching and environmental changes, 
this is not the only perspective from which the role of phase 
variation can be viewed. For example, an alternative sce- 
nario, which does not require a changing environment, was 
recently presented by Ackermann ef aZ. |32|. These authors 
showed that random switching into a "self-sacrificing" phe- 
notypic state can be evolutionarily favoured if the individuals 
in that state have on average greater access to some beneficial 
resource. This idea raises a number of interesting questions 
which we hope to pursue in future research. 

Finally we note that the theoretical framework developed 
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in the present work, although appHed here to the case of detri- 
mental and instantaneous catastrophes, could also be used to 
model environmental changes more generally. For example, 
in the symmetric two state-two environment model consid- 
ered by Thattai and van Oudenaarden |14|, the environment 
flips randomly between two states and these flips are accom- 
panied by a change in fitness from / to 1 — /. This could be 
incorporated in our theoretical framework by setting the /?(/) 
to a constant value and the jump distribution to: 



m(/'I/)-<5(/'-(1-/)) 



(20) 



However, such a choice of fJ,{f\f') would result in fundamen- 
tally different conclusions to those of the present study, be- 
cause the fitness in the Thattai-van Oudenaarden model (and 
in other similar models) is not necessarily decreased when the 
environment changes. In fact, if a large fraction of the cells 
are in the slow-growing state before the environment flips so 
that / < 1/2, then the environmental change will actually in- 
crease the fitness of the population. In contrast, in the present 
work, all catastrophes are detrimental and the advantage of 
switching lies in avoiding the triggering of an environmental 
response. 

The present study suggests a number of avenues for further 
work. First, it would be useful to check the robustness of the 
results to changes in the choice of catastrophe distributions. 
Here we have adopted the power law ([12]) which allows the 
exact solution of the model and generates a broad range of 
catastrophes sizes. Such a distribution could be justified in 
the context of an 'antibiotic' environment, as representing the 
dose-response variability of antimicrobes [331 and variabil- 
ity in the dosage. One could also explore other distributions 
such as exponentially distributed catastrophes or those cen- 
tred about some particular catastrophe fraction /' = af with 
a < 1. It remains to be determined which choice is most 
biologically relevant in different contexts. 

Another point which also deserves investigation in future 
work is the relation between the choice of switching strategy 
and the variability in the population fitness. For example in 
Fig|5]one can see that the different strategies give very differ- 



ent widths for the fitness distribution p{f). In this work we 
defined the optimal strategy as that which gives the maximal 
average growth of the population. However, it might also be 
relevant to include fitness fluctuations in the criteria for opti- 
mality. 

It is also important to consider the case where the environ- 
mental response depends on the absolute size of a particular 
subpopulation. Here, we expect that the population size may 
reach a steady state governed by the balance between growth 
and catastrophes. The total population size could then be max- 
imised either by maximising the growth rate, regardless of 
catastrophes, or by tuning the population composition to avoid 
triggering catastrophes. We thus expect that the two strategies 
identified in the present work will prove to be relevant to a 
variety of models. Moreover, we note that the distinction be- 
tween models based on growth rate and those based on popu- 
lation size may vanish for scenarios with constant population 
size such as chemostat cultures f34]. Equally interesting are 
the prospects for including spatial effects, such as adhesion 
to host surfaces, or transfer between different environmental 
compartments, in the model, and for generalising the model 
to include many different microbial states, in which case the 
same theoretical framework could perhaps be used to describe 
genetic evolution of microbial populations in catastrophic re- 
sponsive environments. 
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Appendix A: Explicit form of p(f) 

Below we give the explicit form for the integral appearing 
in ( 18 1 when /?(/) is given by (|6|: 



df 



m 

v{f) 2A7A/ 



Ice 



(/-/-) 
(/+ - /) 



2A/((/* - /)(r - /-) + + g{f, DgiUJ*)) 
{f-.f-){f*-f-)9if-J*) 



(f -J-) 

s(/_,/*) 



2A/((/ - /*)(/+ - /*) + + g{f, ng{f+J*y 

{f+- f)u+- f*uf+,n 



(/+-/*) 



(Al) 



where 

g{a,h) = ^{a-bf + \^ . (A2) 

From this result, the explicit expression for the fitness distri- 
bution function p{f) can be easily derived. 
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